3D Magnetic Data Inversion using Genetic Algorithms

Table of Contents

• Introduction
• What is 3D Magnetic Data Inversion?
• What is Genetic Algorithm?
• Advantages over Least Square Methods to Find Global Minima
• Applications to Various Cases (Mineral, Oil and Gas, Coal Exploration)

---

Introduction

Magnetic data inversion is a technique used in geophysics to reconstruct the subsurface's magnetic properties by interpreting surface magnetic measurements. Inversion of magnetic data into 3D models helps geophysicists understand subsurface geological structures. Traditional methods for data inversion, like least-square minimization, often get stuck in local minima, making it difficult to find the optimal solution.

To overcome this issue, Genetic Algorithms (GAs) provide a robust, global optimization technique that mimics natural selection to explore a larger solution space. This blog will explore how 3D magnetic data inversion benefits from GAs and their applications in fields such as mineral exploration, oil & gas, and coal mining.

---

What is 3D Magnetic Data Inversion?

3D magnetic data inversion refers to the process of transforming magnetic field data measured at the surface into a 3D model of the subsurface magnetic properties. The aim is to interpret the magnetic anomalies caused by subsurface structures, helping to identify the geometry, depth, and material properties of different geological formations.

Magnetic data inversion is widely used in mineral and hydrocarbon exploration, where magnetic anomalies can indicate the presence of valuable resources. By creating a 3D model of the subsurface, geophysicists can make informed decisions about where to focus exploration efforts.

---

What is Genetic Algorithm?

A Genetic Algorithm (GA) is an evolutionary optimization technique inspired by the process of natural selection. It operates on a population of candidate solutions, evolving them over generations using operators such as selection, crossover (recombination), and mutation.

Key steps in Genetic Algorithms:

1. Initialization: Generate an initial population of solutions.
2. Selection: Choose the best-performing individuals based on a fitness function.
3. Crossover: Combine selected individuals to produce offspring (new solutions).
4. Mutation: Randomly alter some individuals to introduce diversity.
5. Replacement: Replace the old population with the new generation.
6. Termination: Stop when a solution is found or after a set number of generations.

In the context of 3D magnetic data inversion, the fitness function evaluates how well a candidate solution fits the observed magnetic data. Over time, the algorithm converges toward the best-fitting subsurface model.

---

Advantages over Least Square Methods to Find Global Minima

Traditional least-squares methods, such as gradient descent, aim to minimize the difference between observed and predicted magnetic data by adjusting model parameters. However, these methods often converge to local minima due to the non-linear nature of the inversion problem. In contrast, Genetic Algorithms offer several key advantages:

1. Global Search: Genetic Algorithms explore a larger solution space and are less likely to get stuck in local minima.
2. Robustness: GAs do not require the problem to be differentiable, making them ideal for complex, multi-modal landscapes.
3. Diverse Solutions: The population-based approach ensures a diverse set of solutions is explored simultaneously, increasing the chances of finding the global minimum.
4. Adaptability: GAs can handle noisy and incomplete data, which is common in real-world geophysical surveys.

Because of these advantages, GAs are more suitable for solving non-linear inversion problems where traditional methods fail to find the true global solution.

---

Applications to Various Cases (Mineral, Oil and Gas, Coal Exploration)

Genetic Algorithms have been successfully applied to various geophysical exploration problems involving 3D magnetic data inversion, particularly in the following areas:

1. Mineral Exploration

Magnetic surveys are crucial in identifying mineral deposits such as iron ore, nickel, and precious metals. By using GAs, geophysicists can create more accurate 3D models of ore bodies, helping to pinpoint valuable mineral resources and reduce exploration costs.

2. Oil and Gas Exploration

In oil and gas exploration, magnetic data is often used in conjunction with seismic data to identify subsurface structures like faults and sedimentary basins. GA-based inversion methods help generate detailed 3D models that provide insights into potential hydrocarbon reservoirs.

3. Coal Exploration

Magnetic inversion based on aeromagnetic or drone magnetic data can reveal the structure of magnetic dykes present in coalfields, enabling more efficient mining and exploration. GAs are particularly beneficial in coal exploration because they can deal with the complex geology often encountered in coalfields, providing more accurate depth and extent estimates.

---

Conclusion

The application of Genetic Algorithms to 3D magnetic data inversion provides a powerful alternative to traditional least-square methods. Their ability to search globally and handle complex data sets makes them ideal for geophysical exploration in minerals, oil, gas, and coal. By leveraging GAs, geophysicists can create more accurate subsurface models, leading to better decision-making in resource exploration.

Part 1: Creation of Synthetic Magnetic Data to test the inversion algorithm

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import Forward_Model as fm
import time
from tqdm import tqdm

import grav_mag_inv as gm_inv

# Parameters
L = 1000  # Size of the 3D grid (meters)
N = 20   # Number of points in each dimension
cylinder_radius = 200  # Radius of the cylinder (meters)
cylinder_length = 500  # Length of the cylinder along the Z-axis (meters)
cylinder_center = np.array([L/2, L/2])  # Center of the cylinder in the X-Y plane
susceptibility_cylinder = 0.04  # SI
susceptibility_background = 0.001  # SI
susceptibility, X, Y, Z = fm.create_3d_cylinder(cylinder_radius, cylinder_length, susceptibility_cylinder, susceptibility_background, L = L, N = N)

fig = fm.plot_3d_subsurface(susceptibility, X, Y, Z)
fig.show()

# Earth's magnetic field parameters
I = np.radians(49)  # Magnetic inclination in radians
D = np.radians(-10.2)  # Magnetic declination in radians
B_0 = 31000  # Earth's magnetic field strength in nT (assumed)

start_time = time.time()
delta_T_sim = fm.create_3d_forward_model(I, D, B_0, susceptibility, L, N)
elapsed_time = time.time() - start_time
print(f"Time taken: {elapsed_time:.2f} seconds")

Time taken: 0.86 seconds

fm.plot_forward_model(delta_T_sim, X, Y)

Part 2: Inversion to test the algorithm on Synthetic Data

Obsevred Data

t_sim = downsample_observed_data(delta_T_sim, 1.25)

fm.plot_forward_model(t_sim, X_5, Y_5)

# Define F(X) as a function that maps 25000-dimensional X to 2500-dimensional output
def F(X, N):
    #X = X.reshape(N, N, N)
    F_output = fm.create_3d_forward_model(I, D, B_0, X, L, N)
    return F_output

sus_5 = gm_inv.aggregate_to_coarse(X_init, 1.25)

N_5 = int(N/1.25)

# Parameters
population_size = 50  # Number of individuals in the population
epsilon = 0.001  # Small perturbation range

# Generate the initial population
initial_population = gm_inv.generate_initial_population(sus_5, population_size, epsilon)

d_obs = t_sim

initial_population.shape

(50, 16, 16, 16)

N_125 = int(N/1.25)

best_solution, best_error = gm_inv.genetic_algorithm(d_obs, initial_population,population_size=50, N = N_125)

Iterations: 100%|████████████████████████████████████████████████████████████████████| 100/100 [12:57<00:00,  7.77s/it]

best_solution

array([[[0.00148459, 0.00126948, 0.00167638, ..., 0.00090481,
         0.00186301, 0.00088889],
        [0.00135154, 0.00171347, 0.00108971, ..., 0.00147896,
         0.0012531 , 0.00201809],
        [0.0013974 , 0.00220164, 0.00115102, ..., 0.00081912,
         0.00160626, 0.00093439],
        ...,
        [0.00084764, 0.00114032, 0.00113994, ..., 0.00106327,
         0.00071815, 0.00209175],
        [0.00175989, 0.00117551, 0.00131049, ..., 0.00152522,
         0.0012224 , 0.00092644],
        [0.0011119 , 0.00120646, 0.00117466, ..., 0.00218257,
         0.00156473, 0.0024324 ]],

       [[0.00174408, 0.00119189, 0.00058638, ..., 0.00072556,
         0.00241307, 0.00118677],
        [0.00160806, 0.00059292, 0.00184643, ..., 0.00117462,
         0.00215524, 0.00069623],
        [0.00115721, 0.00201371, 0.00131168, ..., 0.00055491,
         0.00198548, 0.00209758],
        ...,
        [0.00102256, 0.00090089, 0.00174717, ..., 0.00151798,
         0.00219178, 0.00192302],
        [0.00047506, 0.00191047, 0.00192699, ..., 0.00177635,
         0.00157967, 0.00140334],
        [0.00243611, 0.0015064 , 0.00175328, ..., 0.00188029,
         0.001434  , 0.00179318]],

       [[0.00228758, 0.00193844, 0.00193965, ..., 0.00147807,
         0.00118974, 0.00142569],
        [0.00235729, 0.00069095, 0.00068143, ..., 0.0011205 ,
         0.00195327, 0.00064193],
        [0.00166034, 0.00148925, 0.0010451 , ..., 0.00215644,
         0.00101173, 0.00153887],
        ...,
        [0.00116507, 0.00094125, 0.00190937, ..., 0.00150599,
         0.00211922, 0.00099319],
        [0.00117636, 0.00253587, 0.0019533 , ..., 0.00115638,
         0.00166693, 0.00071626],
        [0.00140006, 0.00126388, 0.00214869, ..., 0.00178224,
         0.00129882, 0.00231152]],

       ...,

       [[0.00249607, 0.00130331, 0.00084653, ..., 0.00190678,
         0.00138003, 0.00080379],
        [0.00136124, 0.00199141, 0.00251764, ..., 0.00166547,
         0.00086534, 0.00131558],
        [0.00182654, 0.00055228, 0.00185073, ..., 0.00087565,
         0.00070251, 0.00146982],
        ...,
        [0.00218743, 0.00223672, 0.0011305 , ..., 0.00063233,
         0.00108226, 0.00120592],
        [0.000534  , 0.00135732, 0.00133808, ..., 0.00073913,
         0.00203742, 0.00128276],
        [0.00253772, 0.00242554, 0.00185956, ..., 0.00106854,
         0.00129617, 0.00150667]],

       [[0.00205131, 0.00195514, 0.00149625, ..., 0.00074051,
         0.00155021, 0.00199513],
        [0.00057319, 0.00239338, 0.0009035 , ..., 0.00228733,
         0.00081422, 0.00220271],
        [0.00107946, 0.00061654, 0.0023386 , ..., 0.00182296,
         0.00102614, 0.00165934],
        ...,
        [0.00066019, 0.00214824, 0.0018905 , ..., 0.0020003 ,
         0.00077865, 0.00127378],
        [0.00146819, 0.00083949, 0.00147454, ..., 0.00252304,
         0.00154108, 0.00151767],
        [0.00104813, 0.00186259, 0.00083385, ..., 0.0015732 ,
         0.00123245, 0.00085282]],

       [[0.00187474, 0.00166243, 0.00075908, ..., 0.00196238,
         0.00216899, 0.00114306],
        [0.00203486, 0.00151093, 0.00203256, ..., 0.00163873,
         0.0013303 , 0.00238024],
        [0.00152608, 0.00101388, 0.00079974, ..., 0.00216835,
         0.00205742, 0.00197358],
        ...,
        [0.00136234, 0.00216416, 0.00168046, ..., 0.00209954,
         0.00243511, 0.00152112],
        [0.0011474 , 0.00142935, 0.0024032 , ..., 0.0013283 ,
         0.00140646, 0.00179936],
        [0.00083066, 0.00090288, 0.00175912, ..., 0.00153738,
         0.00220227, 0.00180089]]])

Recovered Model

import GravMagPro as gm

fig = fm.plot_3d_subsurface(best_solution, X_5, Y_5, Z_5)
fig.show()

Calculated Data for the recovered Model

start_time = time.time()
delta_T_sim_5 = fm.create_3d_forward_model(I, D, B_0, best_solution, L, int(N/1.25))
elapsed_time = time.time() - start_time
print(f"Time taken: {elapsed_time:.2f} seconds")

Time taken: 0.21 seconds

fm.plot_forward_model(delta_T_sim_5, X_5, Y_5)

Data Misfit

fm.plot_forward_model(delta_T_sim_5-t_sim, X_5, Y_5)

col = gm.get_color()

3D plot of the recovered model with contour surface

fig = gm_inv.plot_3d_subsurface(best_solution, X_5, Y_5, Z_5)
fig.show()

Comparison of True Model and Recovered Model

plt.contourf(X[:,10,:], Z[:,10,:], susceptibility[:,10,:], cmap=col)
plt.colorbar()
plt.title('Original Susceptibility with Cylindrical Anomaly')
plt.xlabel('X (m)')
plt.ylabel('Z (m)')
plt.show()

plt.contourf(X_5[:,10,:], Z_5[:,10,:], best_solution[:,10,:], cmap=col)
plt.colorbar()
plt.title('Recovered Susceptibility with Cylindrical Anomaly')
plt.xlabel('X (m)')
plt.ylabel('Z (m)')
plt.show()